Early stages of radiation damage in graphite and carbon nanostructures: 
A first-principles molecular dynamics study 



o 
o 



in 

(N 



Oleg V. YazyevQ Ivano Tavernelli, Ursula Rothlisberger, and Lothar Helm 
Ecole Polytechmque Federate de Lausanne (EPFL), 
Institute of Chemical Sciences and Engineering, 
CH-1015 Lausanne, Switzerland 
(Dated: February 6, 2008) 

Understanding radiation-induced defect formation in carbon materials is crucial for nuclear tech- 
nology and for the manufacturing of nanostructures with desired properties. Using first principles 
molecular dynamics, we perform a systematic study of the non-equilibrium processes of radiation 
damage in graphite. Our study reveals a rich variety of defect structures (vacancies, interstitials, 
intimate interstitial-vacancy pairs, and in-plane topological defects) with formation energies of 5- 
15 eV. We clarify the mechanisms underlying their creation and find unexpected preferences for 
particular structures. Possibilities of controlled defect-assisted engineering of nanostructures are 
analyzed. In particular, we conclude that the selective creation of two distinct low-energy intimate 
Frenkel pair defects can be achieved by using a 90-110 keV electron beam irradiation. 

PACS numbers: 61.72.Ji, 61.80.Az, 81.05.Uw 
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I. INTRODUCTION 

Radiation resistance of graphite has been one of the 
major concerns of the nuclear industry.— ^ Nowadays, ra- 
diation treatment by high-energy electrons or ions is also 
viewed as a versatile tool for the design of new materials. 
The formation of irradiation-induced defects in graphite- 
like layered carbon nanostructures (multiwalled and bun- 
dled carbon nanotubes, nanoonions, etc.) changes their 
mechanicalS. and electronic properties-^ and may even 
trigger dramatic structural changes. 6,7 However, the 
structure and dynamics of defects in graphite and car- 
bon nanostructures as well as the mechanisms underlying 
their creation and transformation remain elusive. This 
knowledge is crucial for a defect-assisted engineering of 
nanostructures with applications in, e.g., manufacturing 
of nanoelectromechanical systems^ 

Radiation damage of matter is governed by the dis- 
placement of atoms from their equilibrium positions due 
to electronic excitations and direct collisions of high- 
energy particles with the nuclei. In metals and narrow 
band gap semiconductors electronic excitations quench 
instantaneously, leaving collisions with nuclei as the sole 
mechanism responsible for the creation of defects in 
graphite and related carbon materials. 9 If the kinetic 
energy transferred from a high-energy electron or ion 
to the nucleus is higher than the displacement thresh- 
old Td, a carbon atom can leave its initial position to 
form a metastable defect structure on a sub-picosecond 
timescale. Such events are called knock-on displace- 
ments. For highly anisotropic layered carbon materials 
the threshold of the off-plane displacement is T^-«15- 
20 eV— while a creation of defect due to the in-plane 
knock-on collision requires higher transferred energies, 
T^'>30 eV. Possible defects produced by radiation dam- 
age include separated and intimate 1 ^ pairs of intersti- 
tial atoms and vacancies, and in-plane topological de- 
fects involving non-sixmembered rings, e.g. Stone- Wales 



defecti 11 ' 12 The existence of defects in carbon nanostruc- 
tures has been confirmed by direct observations! 13 i 14 

Upon knock-on events a large amount of energy is 
transferred to only a few degrees of freedom. The result- 
ing defect structures formed on a picosecond timescale 
depend on the magnitude and on the direction of the 
transferred momentum and determine the fate of the 
system at longer timescales. Therefore gaining control 
over the early stages of defect formation by tuning the 
irradiation conditions will make the paradigm of the 
defect- assisted engineering feasible. Molecular dynamics 
(MD) simulations performed with empirical potentials 15 
or tight-binding model o 16 d 7 i 18 have been used for the 
studies of radiation damage of various carbon materials. 

In this work, we report a systematic first principles 
study of the early stages of radiation damage of graphite, 
a general model for closely related layered carbon nanos- 
tructures. The paper is organized in the following way: 
In Section [TT] we provide a description of the computa- 
tional methods used in this work. The observed defect 
structures, mechanisms of their formation, and practi- 
cal implications are discussed in Section IIIH Section HVl 
briefly concludes our work. 



II. COMPUTATIONAL METHODS 

By using ab initio molecular dynamics we simulate the 
process of defect formation after the initial transfer of a 
momentum T to one of the carbon atoms in the system. 
The periodic model system consists of a unit cell with 108 
carbon atoms, which contains two graphene sheets with 
stacking ABAB. The dimension of the unit cell in the 
direction perpendicular to the graphene planes was fixed 
to 6.7 A in accordance with the experimental inter-layer 
distance 3.35 A.-— This distance shows only weak varia- 
tion among different layered carbon nanostructures. Our 
computational methodology is based on density func- 
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tional theory (DFT), which lacks a correct description 
of weak van der Waals interactions between graphene 
planes. However, by fixing the unit cell dimension in 
the direction perpendicular to the graphene planes we 
provide a realistic description of layered carbon nanos- 
tructures without any explicit inclusion of van der Waals 
forces. The in-plane distance between two periodic im- 
ages is 12.7 A which is large enough to ensure localization 
of the defect within the unit cell. A coarse sampling of the 
irreducible wedge of the space spanned by the magnitude 
of transferred energy T and the pair of angles <f>E [0° ; 90°] 
and 0e[O°; 60°] (Fig.[[J inset) has been performed. 

The ab initio MD simulations were carried out using 
the CPMD plane wave DFT code 2 ^ and the Perdew, Burke, 
and Ernzerhof exchange-correlation density functional. 21 
A plane wave kinetic energy cutoff of 60 Ry and norm- 
conserving pseudopotentials^ 2 - have been used. The sim- 
ulations were performed within the spin-unrestricted for- 
mulation of DFT starting from an initial guess asym- 
metric with respect to the spin components. Such a 
starting configuration is required in order to ensure a 
broken-symmetry path of bond breaking events^ The 
first 100 fs of each MD simulation were performed us- 
ing the Born-Oppenheimer scheme. The MD timestep 
was set to 0.5 fs. In our simulations we observed that 
during the first 100 fs the transferred kinetic energy was 
well dissipated over the entire system. The initial simu- 
lation was followed by a Car-Parrinello simulation 24 car- 
ried out using a Nose-Hoover thermostat 2 ^, (350 K) un- 
til a stable defect structure was reached (about 1 ps). 
This thermal coupling methodology provides a realistic 
description of the excess kinetic energy dissipation after 
knock-on collisions of reasonably low transferred energies. 
The Car-Parrinello equations of motion were integrated 
with a time step of 0.1 fs using a fictitious electron mass 
of 400 a.u. Finally, the obtained defect structures were 
relaxed by slow annealing of both ionic and electronic 
degrees of freedom. 

The formation energies were evaluated using the 
SIESTA code 2 ^ by relaxing the ionic coordinates and 
the in-plane cell dimensions. The same norm-conserving 
pseudopotentials and density functional as in the plane 
wave calculations together with an optimized double-^ 
plus polarization function (DZP) basis set were used. A 
2x2x2 k-point grid (including the T point) was employed 
in order to obtain accurate defect formation energies* 2 ^ 



III. RESULTS AND DISCUSSION 

A. Off-plane recoils 

The outcomes of our simulations are summarized in 
Figure [Tj (movies of selected MD trajectories are avail- 
able online 28 ). We first discuss the simulation results for 
the off-plane displacements (0e{O°; 30°; 60°}) of carbon 
atoms in inequivalent positions a and /?. The outcomes 
can be divided into four major classes: (i) no defect for- 
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FIG. 1: Polar coordinates representation of the simulation 
outcomes as a function of the parameters T, <j>, and 9. The 
defect structures at each parameter set are given according 
to the nomenclature shown in Figs. [2] and [4] Other labels 
correspond to: 'n' - no defect formation, 'c' - displacement 
cascade, and '-' - simulation not performed. The different 
outcomes corresponding to the off- plane displacements of a//3 
carbon atoms (see left inset) are shown separately. For each 
pair of parameters (4>, 8) the outcomes at different values of T 
are listed according to the values given in the right inset. The 
left inset shows the definition of the parameters determining 
the knock-on displacements. 



mation due to insufficient transferred momentum or due 
to instantaneous recombination of the recoil atom with 
the vacancy ('n'); (ii) separated interstitial- vacancy pairs 
(I-V); (iii) intimate interstitial- vacancy pairs (ilV); (iv) 
displacement cascades (V) in which the recoil atom is 
able to displace other atoms in the lattice. The latter 
case can be viewed as a series of elementary events of 
types (i)-(iii). The simulation of displacement cascades 
is beyond the scope of this study and would require a 
larger unit cell than the one used here. 

The formation of well-separated Frenkel pairs was ob- 
served for atoms in both a and (3 positions at T>25 eV. 
Surprisingly, the interstitial defects were produced only 
in the form of a symmetric "dumbell" structure (I 2 ) 29 i 30 
where the two carbon atoms are symmetrically displaced 
from the graphene plane (Fig. [2j top). Despite the 
highly distorted coordination sphere of these atoms, the 
C-C distance of 1.58 A is close to the one of a typi- 
cal er bond. The core atomic structure is the same as 
for the [l.l.ljpropellane molecule for which a very simi- 
lar C-C bond length (1.60±0.02 A) has been observed 
experimentally^ No single off-plane recoil led to the 
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FIG. 2: Perspective views of the atomic structures of the va- 
cancy (top, left), interstitial (top, right), and intimate Frenkel 
pair (bottom) defects observed in our simulations. The for- 
mation energies are given in parentheses. The values given 
for interstitial defects refer to the formation energies of corre- 
sponding Frenkel pairs. For ilV defects the created vacancy 
is situated in the upper graphene layer. 



"bridge" structure (Ii) 15 i 29 ' 32 with the interstitial atom 
situated between two graphene planes. The formation 
energy of I 2 (£'/=14.3 eV, the value refers to the for- 
mation energy of the corresponding I-V pair) is only 
0.5 eV lower than the one of Ii (Ef = 14.8 eV) where the 
"bridge" interstitial defect is bonded only to the neigh- 
bor atom in the same layer. In this case, a steric re- 
pulsion with the opposite graphene layer contributes to 
the destabilization of the Ii defect. However, bonding to 
the opposite layer leads to more stable shared intersti- 
tial defect structures In the ylid (i5/=14.1 eV) and 
spiro (i?/ = 13.1 eV) configurations, the shared intersti- 
tial atom is additionally bound to one and, respectively, 
to two carbon atoms of the adjacent layer. These struc- 
tures have not been observed in our MD simulations. The 
observed preference for the I2 configuration in graphite 
may have the following origin. In the "dumbell" con- 
figuration the recoil atom is able to transfer its excess 
kinetic energy to the other atoms more efficiently than 
in the case of the "bridge" configuration. At the same 
time, the formation of shared interstitials requires the 
improbable collective motion of a number of atoms in 
the two adjacent graphene layers in the direction of the 
recoil atom. This explains the observed high probabil- 




-0.5 eV 



FIG. 3: The scheme of the self-diffusion process in graphite 
along the c-axis. The diffusing carbon adatom (highlighted in 
the figure) substitutes one of the carbon atoms in the graphite 
layer. The relative energies of local minima and transition 
states are shown. 



ity for the formation of the I2 defect structure in the 
early stages of the radiation-induced defect formation. 
For the isolated graphene sheet, Ii is 0.2 eV more stable 
than I2 due to the absence of the steric repulsion with 
the adjacent graphene layer. In curved graphenic struc- 
tures, like carbon nanotubes, the "bridge" interstitial de- 
fect undergoes further stabilization. Our first principles 
calculations predict that the transition from the defect 
structure I2 to the structure Ii in graphite is character- 
ized by an activation barrier of 0.9 eV. The "dumbell" 
interstitial can also be viewed as a stable intermediate 
of the self-diffusion process in graphite along the c-axis, 
occurring via the substitution of a carbon atom in the 
graphene layer ^ The energy diagram for the diffusion 
process along the c-axis is shown in Figure [3J 

Formation of intimate interstitial-vacancy pairs (ilV) 
requires lower transferred kinetic energies. At T=20 eV 
we observed the formation of two low energy ilV pairs, 
ilVi (_E/=1Q.5 eV) and iIV 2 (E f =ll.O eV) (Fig. [2j bot- 
tom). The displaced atom bridges the defect vacancy 
with two, respectively, three neighbor atoms in the op- 
posite layer, which undergo rehybridization. A fine scan 
of the transferred momentum space indicates a Td value 
of 18 eV for graphite, in agreement with other reported 
values£i2i As a consequence, the use of a particle beam 
energy capable of achieving a maximum kinetic energy 
transfer just above Td will selectively create ilV defects. 
This value for Td would correspond to the maximum ki- 
netic energy transferred by an electron beam of 90 keV. 34 
In the case of carbon nanotubes, Td is expected to be 
lower due to curvature effects This proves the crucial 
role of the ilV defects in the reinforcement of carbon 
nanotube bundles 3 -^ produced by 80 keV electron irra- 
diation. Our results suggest the optimal conditions for 
the modification of mechanic and electronic properties 
of carbon-based layered nanostructures by means of the 
formation of ilV defects. For graphite and closely re- 
lated nanostructures, electron beam acceleration voltages 
of 90-110 kV can be used. Such modifications are non- 
destructive since ilV defects tend to self-recombine with- 
out producing extensive damage of the nanostructure. 14 
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This is also supported by the fact that the barriers for 
ilVi defect recombination™ and for the transformation 
of ilVi into HV2 (0.9 eV in this study) lie below the for- 
mation energies of I-V pairs. 

Our computed formation energy for the previously pro- 
posed ilVi structure^ is in good agreement with the val- 
ues reported in other studies i 10 ' 29 i 35 However, MD simu- 
lations on a longer time scale indicate that the asymmet- 
ric ilVi defect in graphite is not stable against recom- 
bination at 350 K if the shear of neighboring graphite 
layers is allowed. By contrast, the symmetric UV2 de- 
fect is stable throughout our MD simulations. Two 
other intimate Frenkel pairs, 1IV3 (Ef— 11.6 eV) and 1IV5 
(Ef = 12. 2 eV) have been obtained upon off-plane recoils 
caused by larger transferred momenta. In both struc- 
tures the displaced carbon atom is linked to two carbon 
atoms in its host layer and one atom in the neighboring 
layer. It is notable that the formation of ilV defects has 
been observed only upon recoil of the j3 carbon atom. 
We explain this observation by the large probability of 
instantaneous recombination of the a atom recoils due to 
the local arrangement of atoms in the adjacent layer. 



B. In-plane recoils 

The formation of defects after displacement in the 
graphene plane (0=90°) requires higher transferred en- 
ergies T>30 eV. At T=30 eV (0=30°) we observed the 
formation of a Stone- Wales (SW) defect^ which is the 
lowest energy (Ef=4.8 eV) defect in graphite. The mech- 
anism of its formation involves the cyclic permutation of 
three carbon atoms occurring during the first 100 fs after 
the knock-on collision (Fig. [4] top). A much lower activa- 
tion barrier of «10 eV is required when the SW defect is 
formed upon simultaneous in-plane rotation of two neigh- 
boring carbon atomsJ^ However, this mechanism cannot 
be realized upon knock-on collisions because in this case 
the kinetic energy is transferred to a single atom. Ir- 
radiation of graphene-based materials, using an electron 
beam of energy just above 150 keV and oriented along 
the graphene plane, will result in an increase of yield 
of SW defects. This can be used for tuning electronic 
properties of materials^ However, because of the high 
energy transfer required for their formation, SW defects 
will be accompanied by the formation of Frenkel pairs, 
which form upon low energy (T<30 eV) off-plane recoils. 

For T>30 eV two possible general mechanisms of de- 
fect formation have been identified. The first one in- 
volves the formation of strained structures containing 
non-sixmembered rings, which have formation energies 
higher than the formation energy of the SW defect. 
Two such structures, 7-5-4-8 (E f =11.3 eV) and 5-7-7- 
4-4-9 (Ef=12.2 eV) have been observed in our simula- 
tions (Fig. [4l bottom). The second mechanism involves 
the expulsion of one carbon atom from the graphene 
plane shortly after the collision. In this case interstitial- 
vacancy pairs are formed. We observed formation of 
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FIG. 4: Top: The mechanism of formation of a Stone- Wales 
defect upon in-plane knock-on displacement (T=30 eV). The 
carbon atoms involved in the rearrangement are marked 
with letters. Bottom: Atomic structures of 7-5-4-8 and 5- 
7-7-4-4-9 topological defects (arabic numbers indicate non- 
sixmembered rings). 



the "bridge" interstitial defect Ii caused by the expul- 
sion of a carbon atom with low kinetic energy. In ad- 
dition, two new intimate interstitial-vacancy pairs, iIV4 
(E f =12.1 eV) and iIV 6 (Ef =13.0 eV) have been char- 
acterized. In the iIV4 structure the defect is localized 
in the graphene layer where the collision took place. On 
the contrary, the displaced carbon atom in the iIV6 struc- 
ture bridges three atoms of its host layer with one of the 
neighboring layers. The formation energies of all six ilV 
structures found in our simulations lie within a narrow 
interval of 2.5 eV, and they are all below the formation 
energies of separated I-V pairs. These defects should be 
stable at long time scales and at moderate temperatures. 



IV. CONCLUSIONS 

In conclusion, we performed an ab initio molecular dy- 
namics study of radiation-induced defect formation in 
graphite. A variety of different defects, including struc- 
tures which have never been discussed previously, were 
observed in our simulations. The produced defects de- 
pend strongly upon the direction and magnitude of the 
transferred momentum, resulting in the selective forma- 
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tion of certain defect structures. We showed the crucial 
role played by the early stage dynamics in the defect for- 
mation process, and we identified the conditions at which 
selective creation of defects can be achieved. In partic- 
ular, we identified an interval of electron beam energies 
at which only low-energy intimate Frenkel pair defects 
bridging adjacent graphene layers are produced. We also 
conclude that Stone- Wales defects, characterized by the 
lowest formation energy, cannot be produced selectively 
upon irradiation. Our results are of practical importance 
for radiation-assisted manufacturing of carbon materi- 
als and nanostructures with new desired properties and 



functions. 
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